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Summary of Second Year’s Activities 


Below is a brief description of the major activities for this grant during the last year. 
Reprints of two papers (one published and the other in press) describing some of this work 
are attached to this report. 

1. Development of a 200 km—1 Re, () + H + Model. We have modified the basic 
IMTR model so that it deals with a region from 200 km to 1 Re altitude and models 
both the 0+ and H + ions. All of the 0 + -neutral chemistry previously discussed has been 
included as well as the accidentally resonant charge exchange reaction with couples the 0 + 
and H + ions. Currently the electrons are included as a massless neutralizing fluid whose 
density and drift speed match that of the ions and whose temperature profile is given. 

2. Extension of the E x B convection heating study to include centrifugal effects . The 
improved model was used to extend the E x B convection heating study [ Wilson , 1994] 
to higher altitudes so that the consequences of the centrifugal force could be determined. 
Results of this study were presented at the October 1994 Huntsville workshop on Coupling 
of Micro and Mesoscale Processes in Space Plasma Transport and have been submited for 
publication in the proceedings. These results were also presented at the AGU 1994 fall 
meeting in San Francisco. 

5. Study of electron precipitation effects. The model mentioned in l above was used 
to study the effect of energetic electron precipitation. To include the energetic electrons 
we used a two stream electron code developed by Phil Richards. At present the model 
only includes the ionization effects of these electrons and not there effect on thermal 
electron temperature. Preliminary results of this study were presented at the October 
1994 Huntsville workshop and at the 1994 AGU fall meeting in San Francisco. 

4 - Study of wave heating of During 1994 a graduate student (D. Gristina) 

began a study of the consequences of 0 + wave driven heating that occurs at altitudes 
where collisions are important. The waves are produced by precipitating electrons with 
unstable velocity distributions. All the calculations done so far have been local, ignoring 
the effects of transport. Results of this study were presented at the 1994 AGU fall meeting 
in San Francisco. 

5. Polar wind acceleration region study. During 1994 we also made progress in 
data/model comparison study of the polar wind acceleration region. During the summer 
of ? 94 a visiting undergraduate student worked with the RIMS data to find a way to 
sort H+ spin curves by type, principally in terms of the degree of skewing. We now 
know how common the highly skewed spin curves are. Earlier in the year we performed 
a number of simulation runs for comparison with this data. The runs were matched to 
the conditions of the data as much as possible. We found that the model can produced 
H+ velocity distributions which give skewed spin curves, when passed through a RIMS 
instrument simulation filter, but the data often shows more skewing than the current 
model can explain. Results of this study were presented at the 1994 AGU spring meeting 
in Baltimore and at the COSPAR meeting in Hamburg. 
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Plans for the Third Year 

During the coming year we plan to pursue four major areas as described below. 

1. Wave heating induced transport of 0 + . During the coming year D. Gristina 
will use the transport model described in 1 above to study the consequences for 0 + 
transport of wave heating. This will involve modifying the model to include the heating 
algorithm developed for the local calculations discussed in f. above. Parameters used for 
this study will primarily be those appropriate to the TOPAZ III rocket flight. This work 
will constitute Mr. Cristina’s masters thesis. 

2. Proton aurora study. One of the main areas discussed in the original proposal 
which has not been addressed yet is the consequences of energetic ion precipitation for the 
IMTR. We plan, in the next year, to study the effect of proton precipitation, particularly 
the means by which energetic protons could transfer energy to the ionospheric ions. Such 
energy transfer could occur through collisions or plasma instabilities. 

3. Finish study of energetic electron precipitation. The energetic electron precipi- 
tation study will be continued by upgrading the model through the addition of a more 
sophisticated electron description which solves the electron energy equations using expres- 
sion for the electron heat flow. This will allow us to vary the electron temperature with 
time and altitude to see how ion outflow is affect by this aspect of the energetic electron 
precipitation process. 

4 ■ Finish polar cap data-model comparison study. The next step in this study will be 
a statistical analysis of the different types of RIMS spin curves to determine the dependence 
on various geophysical parameters such as altitude, latitude, season, F 10 . 7 , Kp, etc. These 
results will hopefully serve as a guide for additional simulations that might explain the more 
dramatically skewed spin curves. In addition we will study the possibility that a streaming 
instabilities between the H + and 0 + ions in the acceleration region might account for the 
more skewed cases. Recent work done by a collaborator here has suggested a way where 
the microscale effects of the instability can be included in the mesoscale model that is used 
here. 


Papers Published 

Wilson, G. R., Kinetic modeling of 0 + upflows resulting from E x B convection heating 
in the high latitude F region ionosphere, J. Geophys. Res., 99, 17453, 1994. 


Papers Submitted 

* Wilson, G. R., B. L. Giles, and M. 0. Chandler, Collisional to collisionless ion outflow at 
the ionosphere-magnetosphere interface, Adv. Space Res., in press, 1995. 

/Wilson, G. R., The Consequences of Coulomb collisions for space plasma transport, 
in Coupling of Micro and Mesoscale Processes m Space Plasma Transport, AGU 
monograph series, eds. J. L. Horwitz, N. Singh, and J. Burch, 1995. 
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ABSTRACT 


Simple calculations show that ions in the ionosphere are collision dominated while ions in low 
density regions of the magnetosphere are collisionless. For the free outflow of light ions such as H“^~ 
and He“^~, this collisional-collisionless transition occurs in the region where the ions are accelerated 
by the charge separation electric field resulting from the steep gradient in the 0~^” density. In 
this case, specific distortions of the light ion velocity distributions are produced that significantly 
alter ion temperatures and heat flows. The types of nonMaxwellian velocity distributions seen are 
distributions with downward suprathermal tails, inverted bowl shaped distributions, and highly 
skewed distributions with large parallel or perpendicular heat flows. In addition, distributions with 
multiple peaks or one peak and an extended plateau also develop under some conditions. These 
various distributions are the result of relative flow between species, the v'^ dependence of the 
Coulomb collision cross section, and altitude distributed sources of the ions. This paper presents 
results of studies of the evolution of the ion velocity distribution in the transition region done 
with a collisional semikinetic model combining the effects of large scale forces and field-aligned 
transport with the effects of ion-neutral and ion-ion collisions. Comparison of these results with 
data from the RIMS instrument on the DE 1 spacecraft suggest that these distributions have been 
seen and that additional processes, such as downward electric fields, may frequently be present 
on polar cap field lines. 

INTRODUCTION 

It has been known for several decades that the light ions and He+ can freely escape from 
the high latitude ionosphere on field lines open to low pressure regions. Numerous studies of the 
polar wind using a variety of different models (moment based and kinetic) have been performed 
during this time /!/. Some of the basic predictions of the models have been confirmed by direct 
observation / 2 / (and references therein). 

One difficult aspect of polar wind study is understanding how the source of ions in the collision 
dominated ionosphere physically couples to the high speed flow region at high altitudes. This 
requires a knowledge of how ion flow develops while the ions are making a transition between 
collision dominated and collisionless regimes. To study this requires a model which is valid in both 
regimes. The first efforts to accomplish this were by Lemaire / 3 / who combined a hydrodvnamic 
model for the low altitude collisional region with a kinetic model for the high altitude collisionless 
region. A solution was obtained by matching density and velocity moments at the point where 
the two models were joined. Since that time more sophisticated generalized transport models that 
can span the collisional to collisionless transition have been developed /4,5,6/. 

Because the moment based models place restrictions on the form of the velocity distribution, a 
kinetic model is needed to understand how the velocity distribution evolves through the transition 
region. Several recent efforts have been made in this direction / 7 ,8 / inspired bv the Monte, Carlo 
calculations ol Barakat and Lemaire / 9/. Each of these models treats the collision process as a 



random scattering event affecting ion simulation particles as they flow along field lines subject to 
the gravitational, magnetic mirror and ambipolar electric forces. In our model / 7/ the Coulomb 
collision algorithm of Takizuka and Abe / 10/ is used to simulated the effects of Coulomb collisions. 
This technique conserves momentum and energy and gives correct stopping, energy exchange and 
deflection frequencies for the Coulomb collision. 

In this paper we present the results of recent calculations with the semikinetic transition region 
model for the outflow of H^" an d He+ ions in the polar cap in the presence of 0"^ whose 
characteristics are based on DE 1 data. Effort was made in doing these calculations to produce 
results which could be compared with DE 1 RIMS data for the light ions. In addition to 
velocity distributions, the semikinetic model gives altitude profiles of density, drift speed, parallel 
and perpendicular temperature and parallel and perpendicular heat flow. Comparisons between 
density and drift speed from the model are then made with the average of these quantities from 
the light ion data. Additionally, the skew of the RIMS ion spin curves (counts versus spacecraft 
spin phase) is used to compare with synthetic spin curves produced from the model velocity 
distributions. 

RIMS ION DATA IN THE POLAR CAP 

Among the inputs needed for the semikinetic model are the altitude profiles of the 0“^ density 
and temperature as well as the electron temperature. Figure 1 shows RIMS data for 0~^“ density 
versus altitude. The small dots in the figure are individual measurements and the large filled 
circles are mean values of all measurements in a 500 km interval centered on the large circle. The 
individual 0^* density values were obtained from RIMS 0"^ spin curves by the method outlined 
in Chandler et al. /2/. These measurements were taken during the declining phase of the solar 
cycle between December 1981 and December 1983 at invariant latitudes greater than 76°. The 
decrease in the mean values below 1500 km is due to an instrument safety feature which prevented 
it from observing densities much larger than 10 4 cm -3 . The densities are reliable above about 
1700 km. 



The two curves (labeled 4 1’ and ‘2’) in figure 1 are solutions to the O"^ continuity and momentum 
equations assuming quasineutrality with the electrons and a zero 0+ drift speed. For these 
solutions temperature profiles for the electrons and the O - ^ ions were needed. These were obtained 
assuming a constant temperature at low altitude and then an increasing temperature with altitude 
which followed a functional form obtained by assuming that the heat flux was constant and 
proportional to the gradient of the temperature. The constant of proportionality was assumed to 
be proportional to T 5/2 according to Spitzer / 11/. A temperature profile is then determined by its 
values at the two end points and specification of the point, r\ where it goes from being constant to 
being heat flow dominated. For the two curves in figure 1 the temperatures at the lower boundary 
where taken from the IRI-90 model. For curve I the conditions assumed were spring time, 830 





local time. Ap=l8.6 and Fin; = 150 (a large fraction of the data were collected under these 
conditions), lor Curve 2 the conditions assumed were winter time. 830 local time. Ap=I8.6 and 
f ,0.7 = 80. The temperature profiles which correspond to density profile 1 in figure 1 are shown 
in figure 2. The temperatures at the upper boundary and the point r' where adjusted so that the 
density profile matched the mean 0“^ densities from the data. The sum of the electron and 
temperatures at the upper boundary for this case is 15,500 K which is consistent with what is 
needed to match electron densities at higher altitudes assuming that 0+ is the dominant ion to 
higher altitudes over the polar cap / 12/. For density curve 2 the temperatures were assumed to 
be constant with altitude. These density and temperature profiles were used as model inputs. 

In addition to sampling 0“K the RIMS instrument also made measurements of and He~l~. 
Figure 3 shows tour examples of H~^ spin curves at the altitudes and invariant latitudes indicated. 
During the same two year period described above, the RIMS instrument collected about 3500 
1-minute averaged spin curves (like those in figure 3) for both H+ and He+ in the 600 km - 
6378 km altitude range and in the 76° - 90° invariant latitude range. Two bulk parameters which 
can be obtained from these curves are the density and field-aligned flow speed as discussed bv 
Chandler et al, /2/ . In addition, the shape of the spin curve carries information about the 
form of the ion velocity distribution. During the period of data collection the radial head of the 
RIMS instrument operated as an integral instrument. The spin curve then represents a projection 
(or shadow) of the three-dimension velocity distribution as seen from the point in velocity space 
determined by the difference between the spacecraft velocity and the ion bulk velocity. If the 
velocity distribution is symmetric with respect to the magnetic field direction as well as being 
gyrotropic its spin curve will be symmetric from any velocity space vantage point. Skewness of 
the spin curve is an indication of a departure from isotropy in the velocity distribution. Because 
of the geometry of the orbit of DE 1, a negative skewing of the spin curve indicates that the 
velocity distribution is skewed, or elongated, in the downward field-aligned direction. Positive 
skewing indicates the reverse. 
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Fig. 3. Sample DE 1 RIMS H + spin curves 







MODELING ION OUTFLOW IN THE POLAR CAP 


The semikinetic model used here for the transition region calculations is nearly the same as that 
presented in Wilson / 7/ . Both hydrogen and Helium ions, subjected to the ambipolar electric, 
gravitational and magnetic mirror forces, are followed in their one dimensional motion along a 
radial field line. Ion velocities are deflected each time step with random velocity impulses whose 
distribution has a mean of zero and a deviation determined bv the local collision frequency / 10 / . 
Both of these light ions collide with the dominate O^" ion and with themselves while H^-He"^ 
collisions are ignored. This time-dependent model is run until steady-state is reached and then 
time averaged to produce the output. The hydrogen ions are created and removed (mostly in the 
lower portion of the flux tube) by the accidently resonant charge exchange reaction 

H + + 0^H + 0 + (1) 

An additional source of H"^“ is the photoioization of atomic hydrogen. For these reactions the 
model requires altitude profiles of O, H and the neutral temperature (in addition to 0 + ) which 
come from the MSIS -86 model / 13/. The He+ ions are created by photoionization of He and 
removed mainly by reaction with N 2 . These reactions require altitude profiles for He and N 2 
which also come from the MSIS -86 model. 

Figures 4 and 5 show the results for a case (case 1 ) done using the 0^~ density profile from figure 1 
(curve 1 ) and the 0~^ and electron temperature profiles from figure 2. The conditions chosen for 
both the MSIS -86 and IRI-90 models used to generate input for this case were spring time (day 82), 
local time = 830, Ap = 18.6, 80° magnetic latitude, and F 10.7 = 150. The four panels in figure 4 
show the density of all three ion species, the hydrogen and helium field-aligned drift speeds and 
the parallel (‘z 1 ) and perpendicular (‘p’l temperatures and parallel and perpendicular heat flows 
of hydrogen and helium (solid line - H* , dashed line - He^~). The two columns in figure 5 show 
normalized velocity distribution plots for both H^” and He^“ at the indicated altitudes. Each curve 
in the velocity contours is separated by e~ l//2 . Apparent in these results is that 0^~ is the dominant 
ion species at ail altitudes, followed by H+ then He~K The H+ drift speed exceeds the He“t“ drift 
speed at all altitudes and both drift speeds are increasing throughout this altitude range. Both 
light ion species experience primarily perpendicular heating as they collide with the O”^ in the 
1000-2500 km altitude range. In this same region both species have velocity distributions distorted 
in first the upward direction (positive parallel heat flow) and then the downward direction. At 
higher altitudes (above 2500 km) the H"^ ions develop a significant amount of negative skewing 
(negative parallel heat flow) as is apparent in the downward suprathermal tail. This is a result of 
residual H~^ production at high altitude by the charge exchange reaction in equation ( 1 ). The ions 
created at high altitude do not make a significant contribution to the density or drift speed but 
they make a small contribution to the parallel temperature and a large contribution to the parallel 
heat flow. They remain at small field-aligned velocities because they are created at altitudes above 
the region of strong ambipolar electric field acceleration. 

Figures 6 and 7 show results for a case (case 2 ) which uses the second 0“^~ density profile from 
figure 1 . This case is shown as a counter point to that shown in figures 4 and 5 to indicate the 
range of variability in the model results. In the density plots one can see that O^” is no longer the 
dominant ion at all altitudes and that the H~^ and He~^~ densities are nearly the same, especially 
at high altitudes. The light ion drift speeds do not continuously increase with altitude and are 
smaller than the previous case. For both H^~ and He"^ one sees the same basic trend of the parallel 
and perpendicular temperatures at low altitude, but at high altitude the temperature anisotropy 
reverses with the parallel temperature becoming the larger of the two. The hydrogen velocity 
distribution again shows the downward suprathermal tail as seen in the previous case but this tail 
is larger, relative to the core of the distribution. In this case the downward tail is due mainly to 
photoionization of atomic hydrogen at high altitude. This source is now more important relative 
to the charge exchange source because of the smaller amount of 0^” present. The helium ions also 
have a downward tail but in this case the distribution more resembles a Maxwellian core with a 
shoulder on the downward velocity side. These ions are the result of helium ions created at high 
altitude which become trapped in a potential well created on the low altitude side by the strong 
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Fig. 4. Bulk parameters for case 1. 


Fig. 5. Velocity distributions for case 1. 


























Altitude (km) Altitude (km) Altitude (km) 


6000 
5000 
4000 
3000 
2000 
1000 

10° 10 l 10 2 10 3 10 4 10 5 10 6 
Density (cm'O 

6000 
5000 
4000 
3000 
2000 
1000 

0 5 10 15 20 

Drift Speed (km/s) 

6000 
5000 
4000 
3000 
2000 
1000 


0 1000 2000 3000 4000 5000 

Temperature (K) 



Fig. 6. Bulk parameters for case 2. 

























upward ambipoiar electric field (where the electron density gradient is large) and on the topside 
by gravity. At high altitudes the electric field is not strong enough, because of the small density 
gradient and the low electron temperatures, to allow the helium ions created there to overcome 
the gravitational barrier. 

COMPARISON OF MODEL RESULTS WITH TIIE RIMS DATA 

I igure 8 compares the density and drift speed profiles from the two cases discussed in the previous 
section with the mean values of the density and drift speed from the RIMS polar cap data. The 
large filled circles and diamonds give the mean values for each quantity in a 500 km interval 
centered on the point. The solid curves labeled T’ are the corresponding quantities from the 
case shown in figures 4 and 5 while those labeled ; 2’ correspond to figures 6 and 7. It should be 
emphasized here that unlike the case with the 0"^ density in figure 1, there is more scatter of 
the data around the mean values than is present there. These ion data were collected under a 
wide range of solar and geomagnetic conditions during the course of two complete seasonal cycles. 
Never the less, there does appear to be consistent trends in the mean data and there is good 
agreement, at least at low altitude, with the results of case 1. The most important discrepancy 
between the mean data and the model results is the disagreement between the H^~ drift speeds 
at altitudes above 3000 km. A number of other cases were done with different input parameters 
and in each case the model drift speeds were well above the mean of the data. On average, about 
90% of the measured H^ - drift speeds lie well below the model drift speeds above 3000 km. 
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Fig. 8. Comparisons between model results and mean data values. 


Another point of comparison between the data and the model results can be made by comparing 
the shape of the spin curves. To do so we took a model designed to simulate the response of the 
RIMS instrument and adapted it to take as input the model velocity distributions. As additional 
inputs the model requires, the spacecraft velocity in a coordinate system fixed to and aligned with 
the local magnetic field, and the spacecraft potential. The spacecraft velocity is largely constrained 




by the geometry of the orbit of DE L in the polar cap region. Because apogee was well above 1 Re 
altitude the spacecraft was either ascending or descending when it was at altitudes less then 1 Re 
over the polar cap. The spacecraft potential is unknown but estimates based on the densities and 
temperatures from the model suggest that it would range from about -L V at 600 km to about 
2 V at I Re. As an example, figure 9 shows synthetic spin curves for the velocity distribution 
at 4000 km from figure 5. The top row of plots are for an ascending space craft trajectory while 
t he bottom row are for a descending trajectory. The range of assumed spacecraft potentials is 
indicated in each plot. One can see that the downward skewing of the velocity distribution is 
more apparent for an ascending trajectory with a small spacecraft potential. Figure 10 shows 
an ascending trajectory, zero potential spin curve for the H^~ velocity distribution at 4000 km in 
figure 7. More skewing of the spin curve is apparent here because of the large suprathermal tail. 
Some of the spin curves from these two figures resemble the fourth spin curve in figure 3. 
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Fig. 9. Spin curves for H“^ in case I at 4000 km. 
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Fig. 10. Spin curve from case 2. 

Although it is possible to produce spin curves from the model velocity distributions which resemble 
some of the spin curves in the data, many of the observed spin curves have much more skewing 
than can be accounted for by the model results present thus far. What is needed are velocity 
distributions which are much broader with more extended downward tails or shoulders than are 
seen in the velocity distributions in figures 5 or 7. Figure 11 shows two H+ velocity distributions 












at 4000 km for a case where a downward flux of H^“ ions was assumed at the upper boundary and a 
case where an electric potential jump of 1.4 V was assumed at the upper boundary. As can be seen, 
both distributions have more extensive downward tails then in the previous results. Figures 12 and 
14 give the corresponding spin curves for these two distributions in the same format as figure 10. 
These type of distributions can produce spin curves that resemble some of the more common types 
seen in the data such as the second and third curves in figure 3. The presence of a high altitude 
source ot downward H+ flux or a potential barrier has the added benefit of reducing the H^“ drift 
speeds at higher altitudes (above 3000 km) while leaving the lower altitude velocities unchanged. 
It can also raise the density in the model so that it more closely agrees with the mean densities. 
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Fig. 11. H*^ velocity distributions for a) downward flux and b) 1.4 V potential barrier. 
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Fig. 12. Spin curve for 1 1 a). 


Fig. 13. Spin curve for 1 1 b). 


At present these results only suggest the presence of downflowing H*^” or potential barriers over the 
polar cap. An additional mechanisms which might explain the large spin curve skews is the drag 
on the outflowing H*^ ions that would result from waves driven unstable by the relative streaming 
of the H“*~ and 0“^" ions. However, analysis of the model results using the Penrose criterion shows 
that these distributions are stable against the growth of ion acoustic waves. Analysis of their 
stability to the growth of ion cyclotron waves has not yet been done. It has also been suggested 
that time-dependent convection processes may produce results like those seen in figure 11. The 
semikinetic model is time-dependent so that convection or other time-dependent processes could 
be studied. Such tasks will be pursued in the near future. In an effort to better understand what 
the large spin curve skews mean we are currently doing statistical studies of the data to see how 
the spin curve skews relate to location in the polar cap region or to geophysical conditions. 


CONCLUSIONS 


Results from the semikinetic model presented in this paper show that under realistic steady-state 
conditions, the velocity distributions of the and He“^ ions flowing out of the ionosphere in the 







polar rap can deviate significantly from Maxwellian. Often they develop downward skewing in the 
form of a downward suprathermai tail or shoulder on the core distribution. This results in part 
because of the collisions,! drag these ions experience as thev How upward through the nearlv static 
0 + ions and because of small amounts of residual ion production at altitudes above where most of 
the polar wind acceleration occurs. Because of the resemblance between the spin curves produced 
from such model distributions and those seen in the data we believe that such nonMaxwellian 
velocity distributions exist in this region and are explained by the model. However, we have seen 
many spin curves which are more skewed than those produced by the simple model. These larger 
skews indicated that the velocity distributions are more anisotropic than the model indicates. 
To try to produce such distributions and corresponding spin curves, we tried various processes 
including a downward flux of ions at the top boundary and a potential jump at the top boundary. 
These two effects were able to give spin curves which resemble more closely those common in the 
data. Although these two processes are strong candidates for explaining the spin curves with large 
skews, much more work is needed to firmly establish this. 
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ionosphere 
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Abstract. We report here the results of modeling work aimed 
at understanding the development of ionospheric 0 + field-aligned 
upflows that develop in response to high-latitude ExB drift induced 
frictional heating. The model used is a collisional semikinetic 
model which includes ion-neutral resonant charge exchange and 
polarization collisions as well as Coulomb self-collisions. It also 
includes the process of chemical removal of 0 + as well as all of the 
macroscopic forces: ambipolar electric, gravity, magnetic mirror, 
and centripetal. Model results show the development of several 
types of non-Maxwellian velocity distributions including toroids at 
low altitude, distributions with large heat flow in the perpendicular 
component at intermediate altitudes, and distributions with a separate 
upflowing population or upward superthermal tail at high altitudes. 
Whenever the convection electric field increases from a small value 
(< 25 mV/m) to a large value (100-200 mV/m) in 6 min or less large 
upflows develop with parallel drift speeds which peak (below 1000 km) 
at values between 500 m/s and 2 km/s, parallel fluxes which peak 
between 6.0 X 10 8 and 3.2 X 10 9 cm” 2 s” 1 , and parallel per particle heat 
flows which peak between 8.0 X 10“ 9 and 8.0 X 10" 8 ergs cm/s. The 
higher values in these ranges occur for a cooler neutral atmosphere, 
with a larger convection electric field that is turned on quickly. The 
model produces field-aligned 0 + flow speeds that are larger than 
those produced by a 20-moment generalized transport model but 
smaller then those produced by an isotropic hydrodynamic model for 
comparable values of the convection electric field and convection turn 
on times. The model results compare favorably with some topside 
satellite and radar data. 


1. Introduction 

The convection path and speed of plasma in the 
high-latitude ionosphere is highly dependent on pro- 
cesses occurring in the magnetosphere and the solar 
wind. Unlike the ionosphere at low latitudes that gen- 
erally corotates with the Earth, the high-latitude iono- 
sphere can follow complex trajectories at highly vari- 
able speeds. The convection electric field which drives 
this plasma flow via the ExB drift can have values as 
large as 250 mV/m which gives a convection speed 
of about 5 km/s in the lower F region. This pro- 
cess, along with thermospheric winds, can lead to a 
situation where the ions and neutrals flow at large 
speeds relative to each other in a direction perpen- 


Copyright 1994 by the American Geophysical Union. 

Paper number 94JA01214. 

01 48-0227/94/94 JA-0 12 14$05. 00 


dicular to the magnetic field. Since ion and neutral 
thermal velocities in the lower F region are typically 
less than 1 km/s, it is possible, during episodes of 
enhanced convection electric fields, for the ion-neutral 
relative flow to be supersonic. When this happens 
significant perpendicular ion heating can occur as the 
ions and neutrals collide. As has been known for some 
time, with a large enough ExB drift, resonant charge 
exchange collisions between O and O+ can produce 
toroidal ion velocity distributions [Schunk and Walker , 
1972; St-Maunce and Schunk , 1979; Barakat et a/., 
1983; Winkler et a/., 1992]. 

Frictional heating of 0 + ions at low altitude will lead 
to enhanced 0 + upflows. This mechanism can explain 
at least some of the upflows which are seen by radar 
[Yeh and Foster , 1990; Keating et a/., 1990; Wahlund 
et ai y 1992] and low-altitude satellites [Tsunoda et 
aiy 1989; Loranc et ai , 1991; Anderson et a/., 1991]. 
These data sets provide a great deal of information 
concerning the chief characteristics of these upflows, 
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what they are associated with and where they occur. 
Upward velocities in these events (at altitudes of 
1000 km or less) are typically around a few hundred 
meters per second although values above I km/s are 
not uncommon and as high as 3 km/s values have 
been seen. Generally, the velocities increase with 
altitude. Fluxes in these enhanced upflow events 
range from 10 8 to 10 i0 cm -2 s' 1 with typical values 
near 10 9 cm" 2 s' 1 or a few times 10 9 cm" 2 s -1 . 
Sometimes, but not always, the ion temperature is 
elevated in these events. Typical elevated values 
range from 1500 K to about 10,000 K. There is also 
evidence, in some events, for temperature anisotropy 
at low altitude with T±_ > Tjj. In some events the 
electron temperature is elevated, and in some it is not. 
In some events the bottomside electron densities are 
increased, while in others they are reduced. These 
upflows are seen most often in the auroral zone, 
the cusp region and the subauroral region on the 
nightside in so-called “SAID” (subauroral ion drifts) 
events. They are associated with regions of field- 
aligned current, enhanced convection electric fields, 
convection shears, and regions of energetic ion and/or 
electron precipitation. 

There are other possible mechanisms, besides fric- 
tional heating, which could account for these upflows. 
These include current-driven wave heating [Provxdakes 
and Seyler, 1990], energetic electron precipitation (giv- 
ing elevated electron temperatures and ambipolar elec- 
tric fields [Barakat and Schunk , 1983]), ionization of 
energetic neutral oxygen backsplash that results from 
ring current ion precipitation [Torr and Torr, 1979; 
Kozyra et a/., 1982], the auroral density cavity “vac- 
uum cleaner'’ effect [Singh et a/., 1989], and joule heat- 
ing and subsequent expansion of the neutral atmo- 
sphere [Millward et al, 1993]. In many situations, 
no doubt, several of these processes act in concert to 
produce the net effect. In order to determine how much 
of the observed upflow is a result of frictional heating 
or when an upflow event is dominated by frictional 
heating, it is necessary to study this process alone to 
determine its maximum effect and the characteristics 
of the generated upflow. 

The central question in this paper is: How much O* 
upflow will result from ExB driven frictional heating 
in the F region ionosphere and what characteristics 
will this outflow have (i.e., fluxes, drift speeds, temper- 
atures and heat flows, velocity distribution character- 
istics, etc.)? The answer to this question will depend 
on the interplay of a variety of processes including 
ion-neutral collisions, ion-neutral reactions, ion self- 
collisions and the various macroscopic forces such as 
gravity, the mirror force, the ambipolar electric force, 
and the centrifugal force. The collision and reaction 
processes are highly dependent on ion energy, and as 
such, their overall effect will be sensitive to the type 
of ion velocity distribution present. The effect of the 
mirror force is also dependent on the ion velocity dis- 
tribution. In addition to the reactions, collisions, and 
forces, the time history of the frictional heating event 
will also play a role. 

In the last few years there have been a number of 
studies of the F region plasma flows that develop in 


response to iow 1 aiuiude frictional ' heating. Most oi 
these studies have used various types of hydrodynamic 
or generalized transport models. Gombosi and Killeen 
[1987] used a hydrodynamic model, which included 
most of the low-altitude chemical and collisionai pro- 
cesses, to study the effect of a short-duration (about 
5 min) low-altitude frictional heating event. In this 
event the peak, column-integrated heating rate was 
100 ergs cm -2 s' 1 . In response, an 0 + upflow, with a 
flux of about 10 8 cm' 2 s“ l (normalized to 1000 km) 
developed. Above 5000 km the 0+ flow speed ap- 
proached 2 km/s. Cannata et al. [1988], using the 
same model, studied the response of a convecting flux 
tube following realistic convection paths and subject to 
neutral winds derived from the thermospheric global 
circulation model (TGCM) to estimate the heating 
rates used. They concluded that such heating can per- 
turb the ionosphere and that gradients in the heating 
rate may be more important than the absolute heat- 
ing rate in determining the amount of upflow which 
develops. 

In a series of papers [Sellek et a/., 1991; Moffett et 
al 1991; Heelis et a/., 1993] an isotropic hydrody- 
namic model for a closed L = 4 flux tube is used to 
examine the ionospheric response to SAID events for 
different neutral atmospheres (variation in F 10 . 7 ) and 
for various ramp times for the turn on of convection. In 
each case the convection speed was increased linearly 
from 0 to 2 km/s in periods of time ranging from 
15 s to 10 min. Results of these studies show the O+ 
temperature increasing from a value near 1200 K to 
values of 3200 to 4500 K at 300 km altitude. Below 
400 km field-aligned flows are downward but above 
500 km upward flows develop and reach speeds as high 
as 1.9 km/s when the convection speed is turned on 
rapidly. O"*" fluxes at 1000 km altitude range from 
values of about 2 x 10 8 to about 5 x 10 9 cm -2 s _1 
with increasing F 10 . 7 . The maximum 0“*" field-aligned 
drift speed also increases from solar minimum to solar 
maximum. 

Korosmezey et al. [1992] used a more sophisticated 
20-moment generalized transport model to follow a 
flux tube through a 3 km/s SAID event. They 
showed that significant ion temperature anisotropy 
develops in the low-altitude region and that a wave 
like pulse of velocity travels up the flux tube even 
during a period of constant heating. Below 1000 km 
the 0+ parallel drift velocity never exceeded -1-300 m/s 
for this case. When the convection drift speed was 
increased to 4 km/s, the maximum parallel drift (below 
1000 km) reached -1-400 m/s. One possible reason for 
the relatively low parallel drift speeds (compared to the 
isotropic hydrodynamic results) may be the allowance 
for temperature anisotropy in the 20-moment model 
which would lead to a smaller parallel temperature 
increase and a reduction in the pressure gradient force. 
However, the parallel temperature at 200 km does 
reach a value as high or higher then is seen in the 
hydrodynamic results for a convection speed of 3 km/s. 

The work of Loranc and SUMaunce [1993] addresses 
the study of low-altitude frictional heating from a sig- 
nificantly different perspective. They used a collision- 
less kinetic model to study the type of outflows which 
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develop above the exobase (taken at 500 kmi. To 
simulate a frictional heating event, they increased the 
ion temperature at the lower boundary from iOOO to 
5000 K in a period of from 2 to 50 s, held the tempera- 
ture constant for a few minutes and then returned it to 
the lower value. The results of this work showed very 
large parallel drift velocities, between 2 and 3 km/'s 
(below 1000 km altitude), with upward fluxes exceed- 
ing 10 9 cm’ 2 s'" 1 . Frequently, non-Maxwellian veloc- 
ity distributions with two distinct ion populations were 
produced. 

In comparing the isotropic hydrodynamic model, the 
20-moment generalized transport model and the colli- 
sionless kinetic model, it is apparent that large differ- 
ences in the results (field-aligned flow speed) can be 
found depending on whether the model has isotropic 
or anisotropic temperatures, is collisional or collision- 
less or includes kinetic processes (free streaming, mul- 
tiple populations and velocity dispersion). What is 
needed is a model which allows for anisotropic tem- 
peratures and includes both collisions and kinetic pro- 
cesses. Such a model will allow the determination of 
the effect of these various processes on the outflows 
which develop. In this paper we present the results 
of a study of F region ExB frictional heating using 
a collisional kinetic model which includes all of these 
processes. Details of the model are described in the 
next section, with results in the section following. 

2. Description of the Model 

The model used for this work was derived from a 
semikinetic, time-dependent, plasma transport model 
originally developed for the study of the outflow 
of ionospheric plasma at polar latitudes [Wilson et 
a/., 1990] and modified to study the collisional to 
collisionless transition region [ Wilson . 1992]. In this 
study the model was further modified to include ion- 
neutral collisions and reactions. Only the 0+ ion 
is followed, in a simulation region extending from 
30G km to 1100 km. The semikinetic model differs from 
hydrodynamic or generalized transport models in that 
the ions are described through the behavior of a large 
number of simulation particles rather than with a set of 
velocity moments and transport equations. We do not 
follow the full ion cyclotron motion of the simulation 
ions but instead follow the motion of their gyrocenters. 
This approximation is valid if the gyroperiod is small 
compared to any other timescale in the problem. This 
condition is satisfied by all of the results presented in 
this paper. As the simulation gyrocenters move along 
the magnetic field line (assumed vertical), they are 
acted on in a continuous fashion by several large-scale 
macroscopic forces and in a random fashion by various 
microscopic collision and reaction processes. 

The total field-aligned, macroscopic force acting on 
a singly charged O 4 * simulation ion is given by 

dv\\ r 2 £ 2 

F u = m ~dt = + m9 ^ + + 1 - 5_ 3"i2' n) 

In (1). En is the ambipolar electric field along the 
magnetic field direction, g|| is the component of g 
along B. £V is the magnitude of the perpendicular. 


convection electric field at the lower boundary (r 0 ), 
and Bo is the intensity of the magnetic field at r 0 . 
The rest of the symbols in ( 1 ) have their usual meaning 
(i.e.. t. time, rn. ion mass, magnetic moment. V||, 
spatial gradient along the magnetic field; etc.) The 
ambipolar electric field is found by differentiating the 
ambipolar electric potential. The potential is found 
by assuming quasineutrality between the ions and 
electrons and assuming the electrons are isothermal 
and Boltzmann distributed. With this approach and 
an assumed value of the electron temperature ( T e = 
2000 K). knowledge of the ion density is all that is 
needed to find since 

n(0+) = n e = n e0 exp ^ (2) 

We use a value of 2000 K for the uniform electron 
temperature since we wish to concentrate on the effect 
of ion frictional heating on ion upflows. Of the four 
forces listed in (I) the smallest contribution is made 
by the centrifugal force, at least over the altitude range 
considered here. 

Tlie O 4 * ion will undergo a number of different 
types of ion-neutral collisions. With neutral species 
other than atomic oxygen the collision is primarily 
the polarization collision where the ion interacts with 
the electric dipole that it induces in the neutral. 
With atomic oxygen O 4 * ions will experience both 
polarization and resonant charge exchange (RCE) 
collisions. It is the O+-O RCE collisions that are 
primarily responsible for the formation of toroidal 
velocity distributions. For O+-O collisions the cross 
-section is dominated by polarization collisions at low 
energies ( T r < 400 K for colliding Maxwellians) and 
by RCE collisions at high energies [Knof et ai, 1964; 
Banks and Kockarts. 1973]. 

The method used to treat OT-O collisions is es- 
sentially that used by Barakat et ai [1983]. Dur- 
ing each time step, for each simulation ion, a ran- 
dom neutral oxygen velocity is generated from a non- 
drifting Maxwellian velocity distribution with the neu- 
tial exospheric temperature. When E c is nonzero, an 
additional velocity component is added to the neu- 
tral’s velocity to account for the relative convection 
drift between ions and neutrals. (Strictly speaking 
this drift velocity should be added to the ion s velocity 
but because of the interchangabilitv of reference frames 
the collision calculations proceed identically regardless 
of where the drift velocity is added.) The relative 
speed g between the ion and this neutral is then found. 
With this speed and the atomic oxygen density, a ran- 
dom impact parameter h is generated based on Poisson 
statistics 

, / ln( 1 — r) 

~ ]j~nn{0)gdt ^ 

where dt is the time step and r is a random number 
uniformly distributed between 0 and 1. With this same 
value of g , impact parameters for the polarization and 
the RCE collision ( b r and 6 r ) can be found, knowing 
their respective cross sections [Barakat et ai, 1983] 
and the atomic oxygen density. If 6 < b p and b > 6 r 
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only a polarization collision occurs and the resultant 
change in the ion s velocity is calculated assuming the 
collision is elastic and that the ion-neutral relative 
velocity is deflected through a random angle uniformly 
distributed between 0 and tt. If 6 < b r and b > b p an 
RCE collision occurs, but with a probability of In 
that case the ion and neutral velocities are switched. 
If conditions are such that both collisions occur, then 
both are performed with the order randomly selected 
so that the RCE collision occurs first half of the time. 
This collision algorithm was tested by performing 
local calculations identical to the cases considered in 
Barakat et ai [1983]. Our results agreed with theirs 
to within 1%. 

The above approach to 0 + -0 RCE collisions is 
an approximate method. The completely rigorous 
approach to performing these collisions in a particle 
code would be to find, based on the O velocity 
distribution and the given velocity, the distribution 
of the magnitudes of the relative velocity g (i.e. f(g)) 
and then calculate the overall probability of collision 
during a time step dt by 

P - 1 - exp( - dtn(O) J g<r r (g)f(g)dg') (4) 

We assume that f(g) is normalized so that f f{g)dg = 
1. If a collision occurs than an O velocity should be 
randomly selected based on the collision probability 
of the resulting g value and the O+ velocity would 
be changed accordingly. Instead of doing this the 
approach used in the model replaces the integral in the 
exponential in (4) with g<r r (g)n(0)dt where g is found 
at each time step as discussed above. There remains 
the question of how good or bad this approximation 
is. 

If the velocity of the 0+ simulation particle is well 
outside of the region of velocity space where the O 
velocity distribution dominates, then g will be large 
and variations in g that result from selecting different 
O velocities will be small. The approximation of 
neglecting the integral in (4) is a good one since g 
is nearly constant. This situation occurs frequently 
due to the E x B relative drift between the ions and 
neutrals. Even when the ions develop a toroidal 
velocity distribution most of the ions are in a portion 
of velocity space outside of that occupied by the 
neutrals. When the velocity falls within the 

velocity space region occupied by the neutrals g and 
the calculated collision probability will vary widely. 
However, because particles take many time steps to 
cross a cell the total probability for collision for 
this number of time steps will closely approach the 
correct value since the variability in the probabilities 
calculated at each time step will average out upon 
summation. 

The only other ion-neutral collision considered is 
the polarization collision between 0+ and N 2 . The 
procedure employed is identical to that used for the 
polarization collisions between 0 + and O with the 
appropriate change in the polarizability constant found 
in the cross section. Polarization collisions between 
0+ and other neutrals (0 2 , N, H, etc.) were not 
considered since their combined densities where less 


than the combined O and No density by at least 
an order of magnitude for all situations that we 
considered. In order for the ion-neutral collisions to be 
treated properly the time step must be small enough so 
that the fraction of simulation ions experiencing one or 
more collisions in a time step is small. The time step 
used (0.05 s) was small enough so that at 300 km, for 
the cases using the solar maximum atmosphere and 
when the 0 + temperature was highest, < 10% of the 
ions experience one or more collisions in a time step. A 
few tests, performed using half this time step (0.025 s), 
showed no significant differences demonstrating that 
the original time step of 0.05 s is indeed small enough. 
No effort was made to follow multiple ion-neutral 
collisions during a single time step since this happens 
for only a very few particles. For the extreme situation 
just mentioned Poisson statistics have 9.7% of the 
particles experiencing only one collision while 0.3% 
experience two or more. 

In addition to collision processes that can energize 
0+ ions when large ion-neutral relative drifts occur, 
0+ will be lost via reactions with N 2 and O 2 

0+ + N 2 — NO+ + N 
. , (5) 

0 + + Oo — ot + o 

During episodes of frictional heating these reactions 
can lead to a substantial change in the F region iono- 
sphere [Schunk et a/., 1975] since the cross sections for 
these reactions are known to be highly dependent on 
the ion-neutral relative energy [Albritton et ai , 1977]. 
In the model we descritized these cross sections and 
used them to determine whether or not a simulation 
ion was lost during a time step. The procedure was to 
generate a random neutral velocity (for either No or 
Oo) so that the ion-neutral relative speed g could be 
found. (This neutral velocity also has the ExB con- 
vection component added.) This speed is used to find 
the appropriate value of the cross section which is used 
to find the probability of reaction during a time step. 
Losing Poisson statistics and a random number r a deci- 
sion is made about whether or not a reaction occurs. If 
it does, the simulation ion is removed. This procedure 
was tested and found to be working properly since the 
reaction rate as a function of temperature for reacting 
Maxwellians was found to agree with rate equations 
given by Schunk and Nagy [1980]. The effect of the 
ions generated by these reactions was not included in 
the model. These ions will influence the 0 + ions near 
the lower boundary by changing the ambipolar electric 
field and by through Coulomb collisions. Both of these 
effects are expected to be small compared O^-neutral 
collisions. The accidentally resonant charge exchange 
reaction with H (another loss process) was also 
included in the model and handled in the same fashion 
as the reactions with No and 0 2 . 

The effects of small angle Coulomb self-collisions 
are included in the model by use of the method of 
Takizuka and Abe [1977]. If the ion time step is small 
compared to any of the collisional relaxation times, 
this method is equivalent to a solution of the Landau 
form of the Fokker-Planck equation. The small time 
steps needed for the ion-neutral collisions discussed 
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above are adequately small for this purpose. For 
0 + -()^ self-colh sions the procedure is to randomly 
pair G+ simulation ions located in the same spatial 
cell and to compute the deflection of their relative 
velocity U according to that which would result from 
the cumulative effect of many small angle collisions 
experienced over the course of a time step. The 
variable 6 = tani©/ 2 ) (© is the angular deflection of 
F ) is chosen randomly from a Gaussian distribution 
with zero mean and a variance (6 2 ) given by 


(i 2 > = 


e 4 n\nA 

— c it 

2t€ 0 m 2 \U\ 3 


( 6 ) 


where n is the 0+ density in the cell, A is the plasma 
parameter, and m is the oxygen ion mass. The result- 
ing change in the relative velocity U is applied to both 
particles so that between them, momentum and en- 
ergy is conserved. Since these quantities are conserved 
between each colliding pair, they are conserved for the 
entire population. Tests of this procedure show that 
it conserves energy and momentum globally as well as 
relaxing all velocity distributions to the Maxwellian 
over the appropriate timescales. In addition, this pro- 
cedure produces, on average, stopping, deflection and 
energy exchange frequencies as functions of ion veloc- 
ity, appropriate for the Coulomb collision. 

All of the cases that will be discussed in this paper 
were done for dayside (1200 LT) conditions. Since the 
flux tubes were in sunlight, photoionization of O occurs 
throughout each simulation run. This effect is included 
through the introduction of simulation ions in each 
cell at a rate determined by the local atomic oxygen 
density and the ionization rate determined by solar 
conditions (1.7x10“' s “ 1 for solar minimum and 4.6 x 
10 “' s _1 for solar maximum). When new photoions 
are introduced, they are given a randomly generated 
velocity drawn from a Maxwellian distribution that is 
at rest with respect to the neutrals. During the course 


of each run this process only had a small effect since 
each run typically lasted only a few minutes. 

As inputs to the model, a neutral atmosphere with 
altitude profiles of O, H, Oo, No, and the neutral 
temperature is needed. This is obtained from the 
MSIS -86 model of Hedin [1987], Figure 1 shows the 
two neutral atmospheres used for the cases discussed 
in this paper. The MSIS -86 input conditions used were 
the following: local time = 1200, Ap — 60, day = 
183, magnetic latitude = 74°, magnetic longitude = 
180°. For the solar minimum atmosphere a value of 
60 x 10 “ 22 W m “ 2 Hz “ 1 was used for F 10.7 while a 
value of 190 x 10 “ 22 W m “ 2 Hz ” 1 was used for solar 
maximum. 

In doing each case the model was started with the 
0 + ionosphere in diffusive equilibrium with E c = 0 . 
At t = 0 , E c was changed. In some cases it was 
abruptly (ramp time = 0 ) increased to its final value 
(100 or 200 mV/m) and then held constant through 
the rest of the run. In other cases, E c was increased 
linearly to its final value over a period of several 
min (3-6). It is likely that the abrupt E c cases are 
somewhat unphysical, but they serve to illustrate the 
maximum effect. 

In each run the altitude range (300-1100 km) is 
divided into 60 cells with about 200,000 simulation 
ion particles along the whole flux tube. The num- 
ber of simulation ions is maintained near 200,000 by 
renormalization as the number of real ions either in- 
creases or decreases. (Renormalization is accomplished 
by either increasing or decreasing all of the simula- 
tion particle's normalizations simultaneously and then 
either eliminating or adding simulation particles to 
maintain the total content.) The density just below 
the lower boundary is assumed to be the same as that 
of the first cell so that whenever an ion leaves the 
simulation region through the lower boundary a new 
ion is injected upward there. New ions are continu- 
ally injected downward through the upper boundary 
so as to maintain the initial 0"^ density profile and 



Figure 1 . Plots of atomic oxygen, atomic hydrogen, No and Oo for the “solar maximum" 
and “solar minimum ’ atmospheres used in the model. These profiles come from the MSIS -86 
model [Hedin. 1987). A complete listing of the input parameters can be found in the text. 
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Maxwellian velocity distribution. When ions move 
upward across the top boundary, they are removed 
from the simulation. For these reasons, the simulation 
results are only valid to a point slightly beyond the 
time when the bulk of the upflowing ions cross the 
upper boundary. To go much further would require a 
means of handling ions which cross the upper bound- 
ary and then fail back down. 

Model results are presented in the form of altitude 
profiles of bulk parameters (density, drift velocity, flux, 
parallel temperature, perpendicular temperature, and 
parallel heat flow) and velocity space contour plots at 
selected altitudes and times. Bulk parameters were 
found by appropriate summing of particles in a given 
cell at a given time, then averaged for 50 time steps. 
Particles were also summed onto velocity grids (for the 
velocity space plots) at preselected altitudes and time 
averaged for either 50 or 100 time steps. During these 
averaging times the fastest particles seen would travel 
only about two cells if they experience no collisions. 

In doing this type of modeling an important consid- 
eration is the size of the time step compared to the var- 
ious timescales of the problem. Given the densities at 
300 km for the solar maximum atmosphere the shortest 
ion-neutral collision time (for 0+-0 polarization and 
RCE collisions plus 0+-N2 polarization collisions) is 
0.37 s. For the highest 0+ density (1.1 x 10 5 cm" 3 ) 
and lowest O 4 " temperature (880 K) in the solar min- 
imum cases the shortest self-collision time is 0.37 s. 
Likewise, for the solar maximum cases the shortest self- 
collision time is 0.42 s. The shortest chemical reaction 
time occurs when O 4- is highly heated at low altitudes 
and is about 15 s for the solar maximum cases with 
the temperatures which develop in the model results. 
The longest time step used in these simulations is 0.05 s 
which is short compared to the smallest of these times. 
Also, the longest 0+ gyroperiod is 4 ms which is very 
short compared to the various collision times, thus 
justifying the gyrotropic assumption. 

3. Results for both solar minimum and 
solar maximum atmospheres 

Solar minimum cases 

Figure 2 shows bulk parameter time spectrograms 
for a solar minimum case (with time step dt = 
0.05 s) where the convection electric field was abruptly 
increased from 0 to 200 mV/m. The result of this 
change is seen in all of the moments presented. At 
t = 0 the peak O 4- density is about 1.1 x 10 5 cm -3 
at 300 km altitude. Within about 5 min the location 
of the O 4 * peak has moved upward to about 400 km 
and has dropped in value to about 2 x 10 4 cm -3 . This 
change is due primarily to enhanced 0+ removal via 
chemical loss to No and Oo. While this is going on, the 
0+ density at higher altitudes is increasing due to the 
upflow of ions heated above the region where chemical 
loss is significant. In the parallel velocity panel one can 
see a velocity pulse traveling upward which obtains a 
peak speed of about 1.9 km/s near 1000 km after 4 
min. The peak flux is about 3 x 10 9 cm -2 s” 1 , which 
occurs at a lower altitude and an earlier time. 


In the low-altitude ion-neutral collision dominated 
region the 0 + perpendicular temperature rapidly in- 
creases to about 13,000 K as the O 4 " velocity distri- 
bution becomes toroidal. In the same region, ion- 
neutral polarization collisions are responsible for an 
increase of the parallel temperature to about 5000 K. 
Above the ion-neutral collision dominated region 0 + 
self-collisions become more important and the the O 4 " 
velocity distribution relaxes toward isotropy. In this 
region one sees a decline (with altitude) in T± as well 
as an increase in TJ |. At higher altitudes, self-collisions 
decline in importance, and it is possible to develop 
runaway ions whose velocity is large enough for them 
to become collisionless. Evidence of this is seen in the 
large value of the parallel heat flow which is seen to 
develop near the upper boundary after about 3 min. 
Figure 3 shows the O 4 " velocity distribution at the 
upper boundary at several times during this period 
where one can see a nondrifting population and a free 
streaming population. As the number of upflowing 
particles increases (the slower ions reach the upper 
boundary), the free streaming beam merges with the 
nondrifting population. It is at this time when the 
maximum upward bulk flow speed is seen at the upper 
boundary. 

Figure 4 shows a few examples of the velocity 
distributions at the indicated altitudes and times after 
the convection electric field increase. The distribution 
function seen at 40 s and at 950 km is an 880 K 
Maxwellian and is indicative of what is seen at all 
altitudes before the heating begins. Once heating 
begins the velocity distribution at the lower boundary 
(300 km) becomes toroidal very quickly (< 20 s). At 
higher altitudes (~ 400 km) it takes a bit longer as can 
be seen. The toroidal distributions do not extend much 
beyond 400 km in this case since beyond this altitude 
self-collisions become dominant. Above the toroidal 
region the velocity distribution attempts to relax to 
an isotropic form but is not completely successful in 
doing so. Between 500 and 700 km the distributions 
maintain a triangular shape indicating a significant 
amount of upward perpendicular heat flow. Above 
800 km the O 4- densities are low enough, initially, so 
that the newly arriving high speed ions can free stream. 
Later, the free streaming population merges with the 
nondrifting ions because of the arrival of slower ions 
and the large increase in density which causes about a 
fivefold increase in the self-collision frequency. 

In this case one can see three distinct regions which 
develop and change with time. The low-altitude region 
is one where ion-neutral collisions dominate. In this 
region the velocity distribution becomes toroidal and 
the ion perpendicular temperature is about 3 times 
higher than the parallel temperature. This is also a 
region where chemical loss is significant. Above this 
region ion self-collisions dominate and nearly all traces 
of the toroidal distribution are lost. However, because 
of the continued upward flux of ions and energy the 
distribution does not relax to a Maxwellian. In the 
third region two distinct ion populations can develop, 
one streaming relative to the other, for a period of 
time. The location and extent of these regions will 
depend on the neutral atmosphere, the initial ion 
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Figure 2. Altitude versus time spectrograms for bulk parameters from case 1 in Table 1. 
This case used the solar minimum atmosphere from Figure 1 and a convection electric field of 
200 mV/m with a ramp time of 0. 


densities and the strength of the convection electric 
field. 

The development of the ion free streaming region 
depends on the strength of the convection electric field 
and the 0 + density profile. For a fast O+ ion moving 
at a speed u relative to the bulk of the O+ ions its 
mean free path / is given by 


o 2 0 4 

27re“ m~ u * 

e 4 n In A 


( 7 ) 


where m is the ion mass, n is the ion density, and A 


is the plasma parameter. For an ion density of 3 x 
10 4 cm' 3 a relative speed (u) of 1 km/s gives a mean 
free path of only about 1.5 km while a relative speed 
of ~ 5 km/s gives a mean free path of about 1000 km. 
When E r is near 200 mV/m. a significant number of 
ions are heated to speeds near 5 km/s and, with the 0 + 
densities assumed here, become collisionless with re- 
spect to self-collisions. Their increased speed does not 
give them any advantage with respect to ion-neutral 
collisions however. The polarization collision mean 
free path does not change with increasing ion speed, 
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Figure 3. Velocity distributions at 1100 km altitude at the indicated times for case 1 of 
Table 1 (same as case shown in Figure 2). These distributions are normalized so that their 

peak value is one and successive contours represent a decrease in the function by a factor of 
e -i/2 
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Figure 4. Velocity distributions at the indicated altitudes and times for case 1 of Table 1 
(same as case shown in Figure 2). These distributions are normalized so that their peak value 
is one and successive contours represent, a decrease in the function by a factor of e” 1 ^ 2 . 
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while the RCE collision mean free path increases only 
slightly. To free stream, the ions need to be heated to 
speeds near 5 km/'s, and they must be above the region 
where ion-neutral collisions dominate. Free streaming 
is occurring throughout the region above 500 km but 
does not make a significant contribution to the 0 + 
bulk parameters or velocity distributions unless they 
are in a region where the 0“*“ density is relatively low. 
This occurs at altitudes above 850 km before the bulk 
of the slower moving ions reach that altitude. 

A number of other cases were done using the same 
neutral atmosphere, initial 0+ density profile and time 
step. In the first, E c is abruptly increased from 0 
to 100 mV/m and then held constant. In this case, 
things develop more slowly than the 200 mV/m case. 
The peak O+ upward drift speed reaches a value of 
0.95 km/s, while the peak flux is 1.8 x 10 9 cm' 2 s” 1 . 
The maximum perpendicular temperature is about 
3800 K and the maximum parallel temperature is 
about 2500 K. Toroidal velocity distributions do not 
develop (at least above 300 km). There is also little 
effect of free streaming since few ions are heated to 
speeds near 5 km/s. The free streaming velocity 
distributions seen in Figures 3 and 4 do not develop in 
this case. The maximum parallel heat flow seen then is 
about an order of magnitude less than the 200 mV/m 
case discussed above. 

Since the time needed for the effect of the increase 
in E c to propagate through the system is only a few 
minutes, one would expect that when time scales for 
changes in E c are on this order the results would 
change. To see by how much, two cases, in which 
E c was increased linearly, were done for comparison. 
In the first, E c was increased from 0 to 200 mV/m 
in 3 min and then ramped back to zero in 3 min. 
In the second, E c was increased linearly from 0 to 
200 mV/m in 6 min and then held constant. In 
both of these cases the maximum parallel drift was 
reduced from 1.9 to 1.4 km/s and the maximum flux 
was reduced to 2.8 x 10 9 cm” 2 s” 1 for the first case and 


to 2.2 x 10 9 cm - - ^ -1 for the second. The maximum 
parallel heat flow was reduced to 5.2 x 10 -8 ergs cm/s 
for the first case and to 1.5 x 10~ 8 ergs cm/s for the 
second case. 

Because of uncertainties in the value of the 0+- 
O RCE collision cross section [ Salah , 1993] a case 
was done, similar to the first case discussed above, 
where this cross section was doubled in value. As a 
result, the 0 + heating rate increases and the region of 
toroidal distributions expands upward in altitude. The 
maximum drift speed, flux, and parallel heat flow are 
all reduce from those seen in case 1. This case, as well 
as all of the cases discussed above, are summarized in 
Table 1. 

Solar maximum cases 

To determine the effect of changing the neutral 
atmosphere, the cases shown in Table 1 were redone 
with the solar maximum neutral atmosphere shown in 
Figure 1. Figures 5 and 6 summarize, in the same 
way that Figures 2 and 4 do, the results for the 
solar maximum case which corresponds to case 4 in 
Table 1. In this case the convection electric field was 
increased linearly from 0 to 200 mV/m in 6 min and 
held constant after that. This case uses the same time 
step used in case 1 of Table 1 (0.05 s). Because of 
the ramp time things develop more slowly then they 
do for the case of Figure 1, and one does not see the 
development of the maximum drift speed and flux until 
after 6 min. Contrary to what is seen in Figure 1, 
the maximum parallel heat flow develops near the end 
of the simulation, after the maximum drift speed has 
occurred. High perpendicular temperatures occur as 
high as 700 km altitude where the velocity distribution 
becomes toroidal near the end of the simulation. At 
no point in this case do we see distributions where free 
streaming ions comprise a significant fraction of the 
total population. This occurs because the long time 
over which heating is increased provides time for lower- 


Table 1. Solar Minimum Cases 


Case 

t r , min 

E c , mV/m 

Vrnax, km/s 

^ImaX) K 

max » cm 2 s 1 

<3||max> er R s cm /s 

I 

0 

200 

1.9 

13,500 

3.2 x 10 9 

8.0 x 10 -8 

2 

0 

100 

0.95 

3,800 

1.7 x 10 9 

0.8 x 10 -8 

3 

3T 3J 

200 

1.4 

13,000 

2.8 x 10 9 

5.2 x 10 -8 

4 

6 

200 

1.4 

13,500 

<7> 

O 
»— 1 

X 

cp 

ed 

1.5 x 10 -8 

5* 

0 

200 

1.7 

13,700 

2.4 x 10 9 

3.5 x 10 -8 

’"The value a — 

2<r 0 for RCE c 

ollisions 




Table 2. Solar Maximum Cases 




Case 

t r , min 

E c , mV/m 

Fmax, km/s 

flmaxi K 

^|}max> cm s 

Q||max. er K s cm / s 

1 

0 

200 

1.5 

15,000 

2.2 x 10 9 

4.2 x 10~ 8 

2 

0 

100 

0.9 

4.900 

1.3 x 10 9 

0.8 x 10~ 8 

3 

3T 3i 

200 

1.2 

14,400 

2.0 x 10 9 

1.7 x 10 -8 

4 

6 

200 

1.2 

14,400 

1.8 x 10 9 

2.0 x 10 -8 

5* 

0 

200 

1.2 

15,000 

1.3 x 10 9 

2.3 x 10 -8 

6* 

0 

100 

0.75 

4,700 

0.9 x 10 9 

0.9 x 10 -8 


‘The value a = ><r 0 for RCE collisions 
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Figure 5. Altitude versus time spectrograms for bulk parameters from case 4 in Table 2. 
This case used the solar maximum atmosphere from Figure 1 and a convection electric field of 
200 mV/m with a ramp time of 6 min. 


speed ions, heated early, to reach higher altitudes 
before the higher-speed ions, heated later, can. At late 
times in this simulation the O+ velocity distribution 
is toroidal from 300 km to 700 km above which it has 
a predominately triangular shape due mainly to the 
upflow of perpendicular energy. 

Table 2 summarizes several cases done with the 
solar maximum atmosphere shown in Figure I. In 
comparing with the solar minimum cases one can see 
that the maximum perpendicular temperatures are 
slightly higher while the maximum parallel drift speed. 


flux, and heat How are all lower. This difference is 
due to the greater extent of the region where ion- 
neutral collisions dominate. Although the higher 
neutral densities produce greater heating rates, they 
also are responsible for a greater frictional force which 
leads to reduced drift speeds and fluxes. At altitudes 
above 1 100 km one would see, for a solar maximum 
atmosphere, the development of larger parallel drift 
speeds and heat Hows as one more fully enters the 
free streaming region where ion-neutral collisions are 
negligible. 
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Figure 6 . Velocity distributions at the indicated altitudes and times for case 4 of Table 2 
(same as case shown in Figure 5). 


In the various solar minimum and solar maximum 
cases presented above, free streaming ions, ions which 
are upflowing and essentially collisionless with respect 
to self-collisions, are produced whenever O+ ions are 
energized to speeds above 4 km/s. Whether or not 
they make a significant contribution to the total ion 
population (below 1000 km) such that they have an 
impact on bulk parameters such as parallel heat How 
or temperature, depends on the neutral atmosphere 
scale height and on the time history of the convection 
electric field. Those cases for which this happens 
are the solar minimum cases 1, 3, and 5, and the 
solar maximum case 1. The presence of significant 
numbers of free streaming ions occurs, in these model 
results, whenever the parallel per particle heat flow 
(/iii) exceeds about 2.5 x 10~ 8 ergs cm/s or where the 


exceeds 0.5. (For case 1 in Table 1 the maximum value 
of this ration is about 1.2.) Each of the cases where 
this happens are ones where the convection electric 
field increases to its final value fairly rapidly (in less 
than 3 min). It is also more likely to happen when the 
neutral atmosphere is cooler since ion-neutral collisions 
become insignificant at lower altitudes. Also, lower O+ 
densities make it more likely for ions to free stream 
because of the reduced self-collision frequency. Each 
time the ratio in 7 exceeds 0.5 the velocity distribution 
is characterized by a separate upflowing 0 + population 
or an extended superthermal tail. 


4. Discussion 

From the results in this paper we can make the 
following general conclusions about what will happen 
to the topside ionosphere whenever the convection 
electric field increases from a small value ( < 25 mV/ m) 
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>o a large value (100-200 rnV/m) in 0 min or less. 
At altitudes below 1000 km the upward parallel drift 
will peak at values between 500 m/s and 2 km/s, the 
parallel flux peaks at values between 6.0 x L0 8 and 
3.2 x 10 9 cm -2 s“ l , and the parallel per particle 
heat flow peaks at values between 8.0 x 10“ 9 and 

8.0 x 10 ~ 8 ergs cm/s. The higher values in these ranges 
occur for a cooler atmosphere and a larger convection 
electric held that is turned on rapidly. 

Before comparing these results with previous mod- 
eling work it is helpful to understand an issue central 
to this topic. During periods of rapid ExB drift the 
ions are heated primarily in the perpendicular direc- 
tion via RCE collisions with atomic oxygen. In order 
to develop upflows this perpendicular energy must be 
converted to parallel energy. The mirror force is a 
candidate for performing this conversion; however, be- 
cause of the relatively large magnetic field scale length 
it cannot produce strong upflows over the distances in- 
volved here. Another candidate process is the ambipo- 
lar electric field, but in these simulations the electron 
temperature remains uniform and constant, while the 
plasma density gradient decreases so that the electric 
potential drop along the field line actually decreases 
with time. As such, this force makes a decreasing 
contribution to the upflows seen. 0 + -neutral polariza- 
tion collisions and 0+-0+ self-collisions can transfer 
perpendicular energy to parallel energy but, of them- 
selves, cannot impart upward momentum to the ions 
(assuming the neutrals have no vertical drift). They 
can produce elevated parallel temperatures at low al- 
titude which greatly enhances the preexisting parallel 
pressure gradient and which in turn drives the upflows. 

To test the contribution of the two isotropizing col- 
lision processes on the development of upflows, case 1 
of Table 1 was redone, once with both polarization 
and self-collisions turned off and once with only self- 
collisions turned off. 0 + -0 RCE collisions were turned 
on in both cases. In the first of these two cases an 
upsurge of O" 1 " ions did develop; however, its maximum 
parallel velocity was 870 m/s, its maximum flux was 

9.1 x 10 8 cm -2 s“* aud the maximum parallel heat flow 
was 3.3 x 10 -9 ergs cm/s. Because of the lack of self- 
collisions, toroidal free streaming populations began to 
flow up through the existing nondrifting Maxwellian 
population but these free streaming ions only attained 
a parallel flow speed of about 2-3 km/s compared to 
the 4-5 km/s of the free streaming populations seen 
in case 1 of Table i. When only self-collisions were 
turned off. a much stronger upflow developed. In this 
case the maximum parallel velocity was 1.5 km/s, the 
maximum flux was 1.7 x 10 9 cm -2 s" 1 and the max- 
imum parallel heat flow was 6.6 x 10~ 8 ergs cm/s. In 
this case, free streaming populations with drift speeds 
of 4-5 km/s developed which were similar to those 
seen in Figure 3, except that they were broader in 
the perpendicular velocity direction, and they showed 
no tendency to merge with the existing population as 
occurs in the last plot in Figure 3. 

The results of these last few cases clearly illustrate 
the importance of isotropizing collisions for the devel- 
opment of 0 + upflows that result from ExB frictional 
heating. Since these upflows are primarily pressure 


gradient driven Hows the magnitude of the upflow will 
depend on the rate at which the pressure gradient in- 
creases. This increase rate depends on the ion heating 
rate and the efficiency with which polarization and self- 
collisions isotropize the distribution. Polarization col- 
lisions act primarily at low altitude in the region where 
ion-neutral collisions dominate, while self-collisions are 
important mainly above this region. In the strong con- 
vection cases, one might think that the heating of the 
ions in the low altitude ion-neutral collision dominated 
region would leave them energetic enough so that they 
are collisionless to the point where isotropization via 
self-collisions would be inefficient. However, when ions 
are nearly collisionless, there exists an intermediate 
stage where the ions do not exchange much energy 
with each other but their velocity vectors are easily 
deflected. As such, the self-collision process can remain 
fairly efficient at isotropizing the distribution, partic- 
ularly just above the ion-neutral collision dominated 
region where 0+ densities are still high. Kinetic details 
of how these collisions operate become important to 
the overall upflow story. 

Prior to initiation of ion heating the 0+ pressure 
profile is such that the ions are in hydrostatic equilib- 
rium. Once heating begins this equilibrium is upset 
and pressure transport begins to relax the increased 
pressure gradient. Pressure transport is accomplished 
in two ways. One is the transmission of the pressure 
upward in a diffusive process through self-collisions. 
The second method is through the free streaming of 
nearly collisionless energetic ions which introduces the 
effect of velocity dispersion into the results. Both of 
these pressure relaxation processes operate simultane- 
ously in these simulations. 

Differences or similarities of the results in this paper 
and the results from previous work will have much 
to do with the isotropization issue and the ion free 
streaming issue. The model used in the works of Sellefc 
et ai. [1991], Moffett et ai. [1991], and Heehs et ai. 
[1993] assumed an isotropic 0+ temperature so that 
the issue of transfer of perpendicular to parallel energy 
is side stepped since the ions are directly heated in the 
parallel direction. Our results differ significantly from 
those produced by the isotropic hydrodynamic model. 
For cases in the present paper where the maximum 
convection speed is 2 km/s ( E c — 100 mV/m) and 
is turned on instantaneously, the maximum 0 + up- 
flow speed is only about 0.9 km/s compared to values 
twice as large from the hydrodynamic results. This 
difference may be due to the smaller O 9- parallel tem- 
perature which develops in the semikinetic model and 
to differences in t he temperature versus altitude profile 
from the two models. The semikinetic model produces 
upflow velocities as large as 1.9 km/s, but only for 
convection speeds near 4 km/s and for atmospheric 
conditions near solar minimum. Also, in our results we 
see greater upflow speeds and fluxes for solar minimum 
conditions than for solar maximum conditions. This is 
opposite to what is seen in the isotropic hydrodynamic 
results. 

In the paper by Kordsmezey et ai. [1992] cases with 
higher convection electric fields (150 and 200 mV/m) 
were done. The model employed was a 20-moment 
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generalized transport model [Gombosi and Rasmussen. 
1991], which allows for temperature anisotropy and 
parallel and perpendicular heat How. One of the re- 
sults of their calculations is that the maximum upward 
drift speed does not exceed 400 m/s for altitudes be- 
low 1000 km. This conclusion is significantly different 
from ours where the maximum upward drift speed can 
approach 1.5 km/s for similar convection electric held 
strengths in a solar maximum atmosphere. Several 
things may explain this difference. In the Korosmezey 
et al. calculations the convection electric held is turned 
on with a 10-min ramp time compared to our maxi- 
mum ramp time of 6 min. In comparing cases 1 and 4 
in Table 2 one can see that increasing the ramp time 
decreases the maximum drift speed. Another differ- 
ence is that we assumed a uniform electron tempera- 
ture of 2000 K. while they had a steep electron temper- 
ature increase (from the neutral exospheric tempera- 
ture to about 4400 K) over the 200 to 1000 km altitude 
range. (This electron temperature profile is probably 
more appropriate for an open rather than a closed flux 
tube because of the assumed downward electron heat 
flow imposed at the upper boundary.) This will reduce 
the ambipolar electric field and reduce somewhat the 
upflow speeds. In order to make a better comparison 
we performed a case with the same ramp time and 
electron temperature profile. A neutral atmosphere 
at solar maximum. 2400 LT and at a latitude corre- 
sponding to an L shell of L = 5 was also used. In this 
case the parallel drift speed near 1100 km reached a 
value of about 750 m/s by the end of the 10-min ramp 
time and a maximum value of 800 m/s 30 s later. It 
appears that kinetic effects (free streaming ions and 
isotropization) accounts for the remaining differences 
between the results of these two model types. 

In the results of Seliek et ai. [1991], Moffett et 
al. [1991], Heelis et ai [1993] and Korosmezey et 
al. [1992] the 0 + drift speeds which develop at low 
altitudes (< 500 km) are downward with magnitudes 
of 50 to 200 m/s. These downflow's are the result of O' 1- 
chemical loss at low alt itude and the shifting of the F 2 
peak to higher altitude. In our results we see basically 
the same thing: downwards flows of 50-200 m/s below' 
500 km. All of these results contrast sharply with the 
results of Anderson et al. [1991] where large upward 
ion velocities were seen at low altitude in SAID events. 

Above about 500 km, for cases done with the solar 
minimum atmosphere in Figure 1 and a 200 mV/m 
convection electric field, our results compare favorably 
with those of Loranc and St-Maunce [1994], at least as 
far as the drift velocity is concerned. The parallel and 
perpendicular temperature versus altitude profiles are 
however, significantly different; we see more structure 
in the parallel temperature and a larger perpendicular 
temperature gradient. Much of this difference is due 
to self-collisions and the fact that most of energy given 
to the ions in the friction heating process goes into 
perpendicular motion resulting in larger temperature 
anisotropies than are assumed in the Loranc and 
St-Maurice model. For the solar maximum cases 
we see more differences with lower drift speeds at 
high altitudes and greater temperature anisotropies 
at the exobase. Many of the abrupt changes in the 


velocity distribution function seen in their results are 
smoothed out or eliminated in our results because 
of collisions. From our results we concluded that 
above an altitude of about 900 km the Loranc and St- 
Maurice model probably gives resonable results since 
ion-neutral collisions are insignificant there and 0+ 
are low enough so that ions streaming with a speed of 
only 2-3 km/s would be collisionless. 

Our results compare well with some of the obser- 
vations. First, the parallel drift speeds and fluxes 
produced in the model (with convection electric fields 
in the 50-200 mV/m range) fall in the same range as 
those of the thermal ion upwelling events studied by 
Tsunoda et al. [1989] using HILAT data. Second, 
the flow velocities and parallel temperatures seen in 
the 100 mV/m cases compare favorably with what is 
seen in a typical “type 1” event discussed by Wahlund 
et ai [1992]. We also develop higher drift speeds 
(> 1 km/s) as are seen in the “type 2” events, but only 
with convection electric fields greater than 150 mV/m. 
In none of our results do we get parallel drift speeds as 
high as 3 km/s which is seen for upward drifts in some 
radar data [Yek and Foster , 1990] and satellite data 
[Loranc et ai , 1991]. If most of the upflow in those 
data are field-aligned then some other process, such 
as electron precipitation, might need to be invoked 
to account for the higher flow speeds. As previously 
noted, we do not see large parallel drift speeds at low 
altitude as reported by Anderson et al. [1991] 
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